Comparison between the ages for the genders and types

Questions

  • What are the differences between the ages for the different combinations of gender and types?
  • Do we observe the same changes as globally?

Age effect - General Questions

  • What are the differences between the ages?
  • Which genes and pathways are differentially expressed between 8w and 52w, between 52w and 104w, between 8w and 104w? Are they the same? Is there a gradient?
  • Are they different for the two genders?
  • Are they different for the two types?

Loads

Libraries and functions

In [1]:
source("load_libraries.R")

Warning message in is.na(x[[i]]):
“is.na() applied to non-(list or vector) of type 'environment'”Warning message in rsqlite_fetch(res@ptr, n = n):
“Don't need to call dbFetch() for statements, only for queries”
==========================================================================
*
*  Package WGCNA 1.63 loaded.
*
*    Important note: It appears that your system supports multi-threading,
*    but it is not enabled within WGCNA in R. 
*    To allow multi-threading within WGCNA with all available cores, use 
*
*          allowWGCNAThreads()
*
*    within R. Use disableWGCNAThreads() to disable threading if necessary.
*    Alternatively, set the following environment variable on your system:
*
*          ALLOW_WGCNA_THREADS=<number_of_processors>
*
*    for example 
*
*          ALLOW_WGCNA_THREADS=4
*
*    To set the environment variable in linux bash shell, type 
*
*           export ALLOW_WGCNA_THREADS=4
*
*     before running R. Other operating systems or shells will
*     have a similar command to achieve the same aim.
*
==========================================================================


Allowing multi-threading with up to 4 threads.
[1] "preparing gene to GO mapping data..."
[1] "preparing IC data..."
[1] "preparing gene to GO mapping data..."
[1] "preparing IC data..."
[1] "preparing gene to GO mapping data..."
[1] "preparing IC data..."
In [2]:
source("functions.R")

Data

In [3]:
load("../results/dge/gene_length.RData")
load("../results/dge/metadata.RData")
load("../results/dge/norm_counts.RData")
load("../results/dge/dge.RData")
In [4]:
load("../results/dge/dge_net.RData")
load("../results/dge/dge_layout.RData")
load("../results/dge/dge_net_connected_gene_colors.RData")
load("../results/dge/dge_net_pal2.RData")
In [5]:
module_nb = length(unique(connected_gene_colors))
#pal2 = c(pal2, "white", "black")
In [6]:
# Interactions between ages, types and genders
F_SPF_52w_8w = results(dge,contrast= c(0,0,0,0,1,0,0,0,0,0), alpha=0.05, test="Wald")
F_GF_52w_8w = results(dge,contrast= c(0,0,0,0,1,0,0,0,0,1), alpha=0.05, test="Wald")
M_SPF_52w_8w = results(dge,contrast= c(0,0,0,0,1,0,1,0,0,0), alpha=0.05, test="Wald")
M_GF_52w_8w = results(dge,contrast= c(0,0,0,0,1,0,1,0,0,1), alpha=0.05, test="Wald")
F_SPF_104w_8w = results(dge,contrast= c(0,0,0,1,0,0,0,0,0,0), alpha=0.05, test="Wald")
F_GF_104w_8w = results(dge,contrast= c(0,0,0,1,0,0,0,0,1,0), alpha=0.05, test="Wald")
M_SPF_104w_8w = results(dge,contrast= c(0,0,0,1,0,1,0,0,0,0), alpha=0.05, test="Wald")
M_GF_104w_8w = results(dge,contrast= c(0,0,0,1,0,1,0,0,1,0), alpha=0.05, test="Wald")
F_SPF_104w_52w = results(dge,contrast= c(0,0,0,1,-1,0,0,0,0,0), alpha=0.05, test="Wald")
F_GF_104w_52w = results(dge,contrast= c(0,0,0,1,-1,0,0,0,1,-1), alpha=0.05, test="Wald")
M_SPF_104w_52w = results(dge,contrast= c(0,0,0,1,-1,1,-1,0,0,0), alpha=0.05, test="Wald")
M_GF_104w_52w = results(dge,contrast= c(0,0,0,1,-1,1,-1,0,1,-1), alpha=0.05, test="Wald")
In [7]:
to_comp = c("52w VS 8w (F, SPF)","52w VS 8w (F, GF)", "52w VS 8w (M, SPF)", "52w VS 8w (M, GF)",
            "104w VS 52w (F, SPF)", "104w VS 52w (F, GF)", "104w VS 52w (M, SPF)", "104w VS 52w (M, GF)",
            "104w VS 8w (F, SPF)", "104w VS 8w (F, GF)", "104w VS 8w (M, SPF)", "104w VS 8w (M, GF)")
In [10]:
tga_col_order = c(grep("SPF_8w_M_+", colnames(norm_counts), perl=TRUE, value=TRUE),
                  grep("SPF_52w_M_+", colnames(norm_counts), perl=TRUE, value=TRUE),
                  grep("SPF_104w_M_+", colnames(norm_counts), perl=TRUE, value=TRUE),
                  grep("SPF_8w_F_+", colnames(norm_counts), perl=TRUE, value=TRUE),
                  grep("SPF_52w_F_+", colnames(norm_counts), perl=TRUE, value=TRUE),
                  grep("SPF_104w_F_+", colnames(norm_counts), perl=TRUE, value=TRUE),
                  grep("GF_8w_M_+", colnames(norm_counts), perl=TRUE, value=TRUE),
                  grep("GF_52w_M_+", colnames(norm_counts), perl=TRUE, value=TRUE),
                  grep("GF_104w_M_+", colnames(norm_counts), perl=TRUE, value=TRUE),
                  grep("GF_8w_F_+", colnames(norm_counts), perl=TRUE, value=TRUE),
                  grep("GF_52w_F_+", colnames(norm_counts), perl=TRUE, value=TRUE),
                  grep("GF_104w_F_+", colnames(norm_counts), perl=TRUE, value=TRUE))
tga_annot_col = as.data.frame(colData(dge)[, c("age","gender", "type")])
tga_annot_col$age = factor(tga_annot_col$age,c("8w", "52w", "104w"))

gta_col_order = c(grep("SPF_8w_M_+", colnames(norm_counts), perl=TRUE, value=TRUE),
                  grep("SPF_52w_M_+", colnames(norm_counts), perl=TRUE, value=TRUE),
                  grep("SPF_104w_M_+", colnames(norm_counts), perl=TRUE, value=TRUE),
                  grep("GF_8w_M_+", colnames(norm_counts), perl=TRUE, value=TRUE),
                  grep("GF_52w_M_+", colnames(norm_counts), perl=TRUE, value=TRUE),
                  grep("GF_104w_M_+", colnames(norm_counts), perl=TRUE, value=TRUE),
                  grep("SPF_8w_F_+", colnames(norm_counts), perl=TRUE, value=TRUE),
                  grep("SPF_52w_F_+", colnames(norm_counts), perl=TRUE, value=TRUE),
                  grep("SPF_104w_F_+", colnames(norm_counts), perl=TRUE, value=TRUE),
                  grep("GF_8w_F_+", colnames(norm_counts), perl=TRUE, value=TRUE),
                  grep("GF_52w_F_+", colnames(norm_counts), perl=TRUE, value=TRUE),
                  grep("GF_104w_F_+", colnames(norm_counts), perl=TRUE, value=TRUE))
gta_annot_col = as.data.frame(colData(dge)[, c("age","type", "gender")])
gta_annot_col$age = factor(gta_annot_col$age,c("8w", "52w", "104w"))

Differentially expressed genes

In [11]:
age_type_gender_data = list(F_SPF_52w_8w, F_GF_52w_8w, M_SPF_52w_8w, M_GF_52w_8w,
                            F_SPF_104w_52w, F_GF_104w_52w, M_SPF_104w_52w, M_GF_104w_52w,
                            F_SPF_104w_8w, F_GF_104w_8w, M_SPF_104w_8w, M_GF_104w_8w)
names(age_type_gender_data) = to_comp
age_type_gender_deg = extract_diff_expr_genes(age_type_gender_data, "age-effect/age_type_gender/")

Warning message in pcls(G):
“initial point very close to some inequality constraints”Warning message in pcls(G):
“initial point very close to some inequality constraints”Warning message in pcls(G):
“initial point very close to some inequality constraints”Warning message in pcls(G):
“initial point very close to some inequality constraints”Warning message in pcls(G):
“initial point very close to some inequality constraints”Warning message in pcls(G):
“initial point very close to some inequality constraints”Warning message in pcls(G):
“initial point very close to some inequality constraints”Warning message in pcls(G):
“initial point very close to some inequality constraints”Warning message in pcls(G):
“initial point very close to some inequality constraints”Warning message in pcls(G):
“initial point very close to some inequality constraints”Warning message in pcls(G):
“initial point very close to some inequality constraints”Warning message in pcls(G):
“initial point very close to some inequality constraints”Warning message in stack.default(getgo(rownames(l$sign_fc_deg), "mm10", "geneSymbol")):
“non-vector elements will be ignored”Warning message in stack.default(getgo(rownames(as.data.frame(l$deg)), "mm10", "geneSymbol", :
“non-vector elements will be ignored”

Stats

In [13]:
age_type_gender_deg$stat
plot_stat_mat(age_type_gender_deg$stat)
All DEG (Wald padj < 0.05)All over-expressed genes (Wald padj < 0.05 & FC > 0)All under-expressed genes (Wald padj < 0.05 & FC < 0)DEG (Wald padj < 0.05 & abs(FC) >= 1.5)Over-expressed genes (Wald padj < 0.05 & FC >= 1.5)Under-expressed genes (Wald padj < 0.05 & FC <= -1.5)
52w VS 8w (F, SPF)1604 937 667 866 590 276
52w VS 8w (F, GF) 950 601 349 503 369 134
52w VS 8w (M, SPF) 294 210 84 211 177 34
52w VS 8w (M, GF) 189 112 77 148 98 50
104w VS 52w (F, SPF) 280 159 121 178 135 43
104w VS 52w (F, GF) 267 146 121 181 109 72
104w VS 52w (M, SPF)3295177815171300 493 807
104w VS 52w (M, GF)2676135513211366 3551011
104w VS 8w (F, SPF)1502 890 612 911 704 207
104w VS 8w (F, GF)1537 849 688 959 607 352
104w VS 8w (M, SPF)38652077178819971004 993
104w VS 8w (M, GF)2762144013221609 5641045

All DEG (Wald padj < 0.05)

In [14]:
# Differentially expressed genes
upset(as.data.frame(age_type_gender_deg$deg),nsets = 6)

DEG (Wald padj < 0.05 & abs(FC) > 1.5)

In [15]:
upset(as.data.frame(1*(!is.na(age_type_gender_deg$sign_fc_deg))),nsets = 6)

DEG (Wald padj < 0.05 & abs(FC) > 1.5)

Log2FC

In [16]:
fc_annot = data.frame(comp = c(rep("52w VS 8w",4), rep("104w VS 52w",4), rep("104w VS 8w",4)),
                     gender = rep(c(rep("F",2),rep("M",2)),3),
                     type = rep(c("SPF","GF"),6))
rownames(fc_annot) = colnames(age_type_gender_deg$sign_fc_deg)
plot_fc_heatmap(age_type_gender_deg$sign_fc_deg, fc_annot)

Z-score

Column order: type - gender - age

In [18]:
plot_z_score_heatmap(z_scores,
                     rownames(age_type_gender_deg$sign_fc_deg),
                     tga_col_order,
                     tga_annot_col,
                     "All DE genes in comparison between the ages for the genders and types",
                     tga_col_order)

comps = list(
    "52w VS 8w (F, SPF)" = c(grep("SPF_8w_F_+", colnames(norm_counts), perl=TRUE, value=TRUE),
                  grep("SPF_52w_F_+", colnames(norm_counts), perl=TRUE, value=TRUE)),
    "104w VS 52w (F, SPF)" = c(grep("SPF_52w_F_+", colnames(norm_counts), perl=TRUE, value=TRUE),
                  grep("SPF_104w_F_+", colnames(norm_counts), perl=TRUE, value=TRUE)),
    "104w VS 8w (F, SPF)" = c(grep("SPF_8w_F_+", colnames(norm_counts), perl=TRUE, value=TRUE),
                  grep("SPF_104w_F_+", colnames(norm_counts), perl=TRUE, value=TRUE)),
    "52w VS 8w (M, SPF)" = c(grep("SPF_8w_M_+", colnames(norm_counts), perl=TRUE, value=TRUE),
                  grep("SPF_52w_M_+", colnames(norm_counts), perl=TRUE, value=TRUE)),
    "104w VS 52w (M, SPF)" = c(grep("SPF_52w_M_+", colnames(norm_counts), perl=TRUE, value=TRUE),
                  grep("SPF_104w_M_+", colnames(norm_counts), perl=TRUE, value=TRUE)),
    "104w VS 8w (M, SPF)" = c(grep("SPF_8w_M_+", colnames(norm_counts), perl=TRUE, value=TRUE),
                  grep("SPF_104w_M_+", colnames(norm_counts), perl=TRUE, value=TRUE)),    
    "52w VS 8w (F, GF)" = c(grep("GF_8w_F_+", colnames(norm_counts), perl=TRUE, value=TRUE),
                  grep("GF_52w_F_+", colnames(norm_counts), perl=TRUE, value=TRUE)),
    "104w VS 52w (F, GF)" = c(grep("GF_52w_F_+", colnames(norm_counts), perl=TRUE, value=TRUE),
                  grep("GF_104w_F_+", colnames(norm_counts), perl=TRUE, value=TRUE)),
    "104w VS 8w (F, GF)" = c(grep("GF_8w_F_+", colnames(norm_counts), perl=TRUE, value=TRUE),
                  grep("GF_104w_F_+", colnames(norm_counts), perl=TRUE, value=TRUE)),
    "52w VS 8w (M, GF)" = c(grep("GF_8w_M_+", colnames(norm_counts), perl=TRUE, value=TRUE),
                  grep("GF_52w_M_+", colnames(norm_counts), perl=TRUE, value=TRUE)),
    "104w VS 52w (M, GF)" = c(grep("GF_52w_M_+", colnames(norm_counts), perl=TRUE, value=TRUE),
                  grep("GF_104w_M_+", colnames(norm_counts), perl=TRUE, value=TRUE)),
    "104w VS 8w (M, GF)" = c(grep("GF_8w_M_+", colnames(norm_counts), perl=TRUE, value=TRUE),
                  grep("GF_104w_M_+", colnames(norm_counts), perl=TRUE, value=TRUE))
)

for(comp in names(comps)){
    plot_z_score_heatmap(z_scores,
                     rownames(age_type_gender_deg$sign_fc_deg)[!is.na(age_type_gender_deg$sign_fc_deg[,comp])],
                     tga_col_order,
                     tga_annot_col,
                     paste("DE genes in", comp),
                     comps[[comp]])
}

Column order: gender - type - age

In [19]:
plot_z_score_heatmap(z_scores,
                     rownames(age_type_gender_deg$sign_fc_deg),
                     gta_col_order,
                     gta_annot_col,
                     "All DE genes in comparison between the ages for the genders and types",
                     gta_col_order)

comps = list(
    "52w VS 8w (F, SPF)" = c(grep("SPF_8w_F_+", colnames(norm_counts), perl=TRUE, value=TRUE),
                  grep("SPF_52w_F_+", colnames(norm_counts), perl=TRUE, value=TRUE)),
    "52w VS 8w (F, GF)" = c(grep("GF_8w_F_+", colnames(norm_counts), perl=TRUE, value=TRUE),
                  grep("GF_52w_F_+", colnames(norm_counts), perl=TRUE, value=TRUE)),
    "104w VS 52w (F, SPF)" = c(grep("SPF_52w_F_+", colnames(norm_counts), perl=TRUE, value=TRUE),
                  grep("SPF_104w_F_+", colnames(norm_counts), perl=TRUE, value=TRUE)),
    "104w VS 52w (F, GF)" = c(grep("GF_52w_F_+", colnames(norm_counts), perl=TRUE, value=TRUE),
                  grep("GF_104w_F_+", colnames(norm_counts), perl=TRUE, value=TRUE)),
    "104w VS 8w (F, SPF)" = c(grep("SPF_8w_F_+", colnames(norm_counts), perl=TRUE, value=TRUE),
                  grep("SPF_104w_F_+", colnames(norm_counts), perl=TRUE, value=TRUE)),
    "104w VS 8w (F, GF)" = c(grep("GF_8w_F_+", colnames(norm_counts), perl=TRUE, value=TRUE),
                  grep("GF_104w_F_+", colnames(norm_counts), perl=TRUE, value=TRUE)),
    "52w VS 8w (M, SPF)" = c(grep("SPF_8w_M_+", colnames(norm_counts), perl=TRUE, value=TRUE),
                  grep("SPF_52w_M_+", colnames(norm_counts), perl=TRUE, value=TRUE)),
    "52w VS 8w (M, GF)" = c(grep("GF_8w_M_+", colnames(norm_counts), perl=TRUE, value=TRUE),
                  grep("GF_52w_M_+", colnames(norm_counts), perl=TRUE, value=TRUE)),
    "104w VS 52w (M, SPF)" = c(grep("SPF_52w_M_+", colnames(norm_counts), perl=TRUE, value=TRUE),
                  grep("SPF_104w_M_+", colnames(norm_counts), perl=TRUE, value=TRUE)),
    "104w VS 52w (M, GF)" = c(grep("GF_52w_M_+", colnames(norm_counts), perl=TRUE, value=TRUE),
                  grep("GF_104w_M_+", colnames(norm_counts), perl=TRUE, value=TRUE)),
    "104w VS 8w (M, SPF)" = c(grep("SPF_8w_M_+", colnames(norm_counts), perl=TRUE, value=TRUE),
                  grep("SPF_104w_M_+", colnames(norm_counts), perl=TRUE, value=TRUE)),
    "104w VS 8w (M, GF)" = c(grep("GF_8w_M_+", colnames(norm_counts), perl=TRUE, value=TRUE),
                  grep("GF_104w_M_+", colnames(norm_counts), perl=TRUE, value=TRUE))
)

for(comp in names(comps)){
    plot_z_score_heatmap(z_scores,
                     rownames(age_type_gender_deg$sign_fc_deg)[!is.na(age_type_gender_deg$sign_fc_deg[,comp])],
                     gta_col_order,
                     gta_annot_col,
                     paste("DE genes in", comp),
                     comps[[comp]])
}

Co-expression (WGCNA)

DEG into gene co-expression network

  • White: up-regulated
  • Black: down-regulated

52w VS 8w

52w VS 8w M F
SPF
GF
In [20]:
#par(mfrow=c(2,2),mar=c(0,0,0,0))
#col_52w_vs_8w_F_GF = get_deg_colors(age_type_gender_deg, "52w VS 8w (M, SPF)", connected_gene_colors, module_nb)
#plot_net_with_layout(net, col_52w_vs_8w_F_GF, pal2, layout, add_legend=F)
#col_52w_vs_8w_M_GF = get_deg_colors(age_type_gender_deg, "52w VS 8w (F, SPF)", connected_gene_colors, module_nb)
#plot_net_with_layout(net, col_52w_vs_8w_M_GF, pal2, layout, add_legend=F)
#col_52w_vs_8w_F_SPF = get_deg_colors(age_type_gender_deg, "52w VS 8w (M, GF)", connected_gene_colors, module_nb)
#plot_net_with_layout(net, col_52w_vs_8w_F_SPF, pal2, layout, add_legend=F)
#col_52w_vs_8w_M_SPF = get_deg_colors(age_type_gender_deg, "52w VS 8w (F, GF)", connected_gene_colors, module_nb)
#plot_net_with_layout(net, col_52w_vs_8w_M_SPF, pal2, layout, add_legend=F)

104w VS 52w

104w VS 52w M F
SPF
GF
In [21]:
#par(mfrow=c(2,2),mar=c(0,0,0,0))
#col_104w_vs_52w_F_GF = get_deg_colors(age_type_gender_deg, "104w VS 52w (M, SPF)", connected_gene_colors, module_nb)
#plot_net_with_layout(net, col_104w_vs_52w_F_GF, pal2, layout, add_legend=F)
#col_104w_vs_52w_M_GF = get_deg_colors(age_type_gender_deg, "104w VS 52w (F, SPF)", connected_gene_colors, module_nb)
#plot_net_with_layout(net, col_104w_vs_52w_M_GF, pal2, layout, add_legend=F)
#col_104w_vs_52w_F_SPF = get_deg_colors(age_type_gender_deg, "104w VS 52w (M, GF)", connected_gene_colors, module_nb)
#plot_net_with_layout(net, col_104w_vs_52w_F_SPF, pal2, layout, add_legend=F)
#col_104w_vs_52w_M_SPF = get_deg_colors(age_type_gender_deg, "104w VS 52w (F, GF)", connected_gene_colors, module_nb)
#plot_net_with_layout(net, col_104w_vs_52w_M_SPF, pal2, layout, add_legend=F)

Z-score in modules

In [22]:
mod_pal = pal2
names(mod_pal) = paste("ME", names(pal2), sep='')
names(mod_pal) = replace(names(mod_pal), which(names(mod_pal) == 'ME0'), "No module")
annot_colors = list(
    module = mod_pal
)

Column order: type - gender - age

In [23]:
plot_z_score_heatmap_with_modules(z_scores,
                                  rownames(z_scores),
                                  tga_col_order,
                                  tga_annot_col,
                                  "All genes")

comps = c("52w VS 8w (F, SPF)","104w VS 52w (F, SPF)","104w VS 8w (F, SPF)",
          "52w VS 8w (M, SPF)","104w VS 52w (M, SPF)", "104w VS 8w (M, SPF)",    
          "52w VS 8w (F, GF)", "104w VS 52w (F, GF)", "104w VS 8w (F, GF)",
          "52w VS 8w (M, GF)", "104w VS 52w (M, GF)", "104w VS 8w (M, GF)")
for(comp in comps){
    plot_z_score_heatmap_with_modules(z_scores,
                                     rownames(age_type_gender_deg$sign_fc_deg)[!is.na(age_type_gender_deg$sign_fc_deg[,comp])],
                                     tga_col_order,
                                     tga_annot_col,
                                     paste("DE genes in", comp))
}

Column order: gender - type - age

In [24]:
plot_z_score_heatmap_with_modules(z_scores,
                                  rownames(z_scores),
                                  gta_col_order,
                                  gta_annot_col,
                                  "All genes")
comps = c("52w VS 8w (F, SPF)","52w VS 8w (F, GF)",
          "104w VS 52w (F, SPF)","104w VS 52w (F, GF)",
          "104w VS 8w (F, SPF)","104w VS 8w (F, GF)",
          "52w VS 8w (M, SPF)","52w VS 8w (M, GF)",
          "104w VS 52w (M, SPF)","104w VS 52w (M, GF)",
          "104w VS 8w (M, SPF)","104w VS 8w (M, GF)")
for(comp in comps){
    plot_z_score_heatmap_with_modules(z_scores,
                                     rownames(age_type_gender_deg$sign_fc_deg)[!is.na(age_type_gender_deg$sign_fc_deg[,comp])],
                                     gta_col_order,
                                     gta_annot_col,
                                     paste("DE genes in", comp))
}

Genes in modules

In [28]:
comps = c("52w VS 8w (F, SPF)","52w VS 8w (F, GF)",
          "104w VS 52w (F, SPF)","104w VS 52w (F, GF)",
          "104w VS 8w (F, SPF)","104w VS 8w (F, GF)",
          "52w VS 8w (M, SPF)","52w VS 8w (M, GF)",
          "104w VS 52w (M, SPF)","104w VS 52w (M, GF)",
          "104w VS 8w (M, SPF)","104w VS 8w (M, GF)")
for(comp in comps){
    plot_top_deg_in_modules(age_type_gender_deg$sign_fc_deg, comp, connected_gene_colors)
}
options(repr.plot.width=7, repr.plot.height=7)

GO analysis

In [29]:
full_go_desc = age_type_gender_deg$GO_wall[[1]][,"term"]
names(full_go_desc) = age_type_gender_deg$GO_wall[[1]][,"category"]
comp = colnames(age_type_gender_deg$over_represented_GO)
comp = comp[4:length(comp)]

Biological process

Dot-plot with the most over-represented BP GO (20 most significant p-values for the different comparison)

In [30]:
plot_top_go(age_type_gender_deg,
            "BP",
            40)

Using term, id as id variables
Using term, id as id variables

Network based on description similarity

In [32]:
BP_network = create_GO_network(age_type_gender_deg, "BP", BP_GO)

52w VS 8w

52w VS 8w M F
SPF
GF
In [33]:
par(mfrow=c(2,2),mar=c(0,0,0,0))
plot_GO_networks(BP_network, "52w VS 8w (M, SPF)", full_go_desc, plot_interactive = FALSE)
plot_GO_networks(BP_network, "52w VS 8w (F, SPF)", full_go_desc, plot_interactive = FALSE)
plot_GO_networks(BP_network, "52w VS 8w (M, GF)", full_go_desc, plot_interactive = FALSE)
plot_GO_networks(BP_network, "52w VS 8w (F, GF)", full_go_desc, plot_interactive = FALSE)

52w VS 8w (M, SPF)

In [34]:
col = get_GO_network_col(BP_network, "52w VS 8w (M, SPF)")
dotRes = getAmigoTree(goIDs=names(col),
                      color=col,
                      filename="../results/dge/age-effect/age_type_gender/go/52w_VS_8w_M_SPF",
                      picType="png",
                      saveResult=TRUE)

GO Tree at "../results/dge/age-effect/age_type_gender/go/52w_VS_8w_M_SPF.png"

52w VS 8w (F, SPF)

In [35]:
col = get_GO_network_col(BP_network, "52w VS 8w (F, SPF)")
dotRes = getAmigoTree(goIDs=names(col),
                      color=col,
                      filename="../results/dge/age-effect/age_type_gender/go/52w_VS_8w_F_SPF",
                      picType="png",
                      saveResult=TRUE)

GO Tree at "../results/dge/age-effect/age_type_gender/go/52w_VS_8w_F_SPF.png"

52w VS 8w (M, GF)

In [36]:
col = get_GO_network_col(BP_network, "52w VS 8w (M, GF)")
dotRes = getAmigoTree(goIDs=names(col),
                      color=col,
                      filename="../results/dge/age-effect/age_type_gender/go/52w_VS_8w_M_GF",
                      picType="png",
                      saveResult=TRUE)

GO Tree at "../results/dge/age-effect/age_type_gender/go/52w_VS_8w_M_GF.png"

52w VS 8w (F, GF)

In [37]:
col = get_GO_network_col(BP_network, "52w VS 8w (F, GF)")
dotRes = getAmigoTree(goIDs=names(col),
                      color=col,
                      filename="../results/dge/age-effect/age_type_gender/go/52w_VS_8w_F_GF",
                      picType="png",
                      saveResult=TRUE)

GO Tree at "../results/dge/age-effect/age_type_gender/go/52w_VS_8w_F_GF.png"

104w VS 52w

104w VS 52w M F
SPF
GF
In [38]:
par(mfrow=c(2,2),mar=c(0,0,0,0))
plot_GO_networks(BP_network, "104w VS 52w (M, SPF)", full_go_desc, plot_interactive = FALSE)
plot_GO_networks(BP_network, "104w VS 52w (F, SPF)", full_go_desc, plot_interactive = FALSE)
plot_GO_networks(BP_network, "104w VS 52w (M, GF)", full_go_desc, plot_interactive = FALSE)
plot_GO_networks(BP_network, "104w VS 52w (F, GF)", full_go_desc, plot_interactive = FALSE)

104w VS 52w (M, SPF)

In [39]:
col = get_GO_network_col(BP_network, "104w VS 52w (M, SPF)")
dotRes = getAmigoTree(goIDs=names(col),
                      color=col,
                      filename="../results/dge/age-effect/age_type_gender/go/104w_VS_52w_M_SPF",
                      picType="png",
                      saveResult=TRUE)

GO Tree at "../results/dge/age-effect/age_type_gender/go/104w_VS_52w_M_SPF.png"

104w VS 52w (F, SPF)

In [40]:
col = get_GO_network_col(BP_network, "104w VS 52w (F, SPF)")
dotRes = getAmigoTree(goIDs=names(col),
                      color=col,
                      filename="../results/dge/age-effect/age_type_gender/go/104w_VS_52w_F_SPF",
                      picType="png",
                      saveResult=TRUE)

GO Tree at "../results/dge/age-effect/age_type_gender/go/104w_VS_52w_F_SPF.png"

104w VS 52w (M, GF)

In [41]:
col = get_GO_network_col(BP_network, "104w VS 52w (M, GF)")
dotRes = getAmigoTree(goIDs=names(col),
                      color=col,
                      filename="../results/dge/age-effect/age_type_gender/go/104w_VS_52w_M_GF",
                      picType="png",
                      saveResult=TRUE)

GO Tree at "../results/dge/age-effect/age_type_gender/go/104w_VS_52w_M_GF.png"

104w VS 52w (F, GF)

In [42]:
col = get_GO_network_col(BP_network, "104w VS 52w (F, GF)")
dotRes = getAmigoTree(goIDs=names(col),
                      color=col,
                      filename="../results/dge/age-effect/age_type_gender/go/104w_VS_52w_F_GF",
                      picType="png",
                      saveResult=TRUE)

GO Tree at "../results/dge/age-effect/age_type_gender/go/104w_VS_52w_F_GF.png"

104w VS 8w

104w VS 8w M F
SPF
GF
In [43]:
par(mfrow=c(2,2),mar=c(0,0,0,0))
plot_GO_networks(BP_network, "104w VS 8w (M, SPF)", full_go_desc, plot_interactive = FALSE)
plot_GO_networks(BP_network, "104w VS 8w (F, SPF)", full_go_desc, plot_interactive = FALSE)
plot_GO_networks(BP_network, "104w VS 8w (M, GF)", full_go_desc, plot_interactive = FALSE)
plot_GO_networks(BP_network, "104w VS 8w (F, GF)", full_go_desc, plot_interactive = FALSE)

104w VS 8w (M, SPF)

In [44]:
col = get_GO_network_col(BP_network, "104w VS 8w (M, SPF)")
dotRes = getAmigoTree(goIDs=names(col),
                      color=col,
                      filename="../results/dge/age-effect/age_type_gender/go/104w_VS_8w_M_SPF",
                      picType="png",
                      saveResult=TRUE)

GO Tree at "../results/dge/age-effect/age_type_gender/go/104w_VS_52w_M_SPF.png"

104w VS 8w (F, SPF)

In [45]:
col = get_GO_network_col(BP_network, "104w VS 8w (F, SPF)")
dotRes = getAmigoTree(goIDs=names(col),
                      color=col,
                      filename="../results/dge/age-effect/age_type_gender/go/104w_VS_8w_F_SPF",
                      picType="png",
                      saveResult=TRUE)

GO Tree at "../results/dge/age-effect/age_type_gender/go/104w_VS_52w_F_SPF.png"

104w VS 8w (M, GF)

In [46]:
col = get_GO_network_col(BP_network, "104w VS 8w (M, GF)")
dotRes = getAmigoTree(goIDs=names(col),
                      color=col,
                      filename="../results/dge/age-effect/age_type_gender/go/104w_VS_8w_M_GF",
                      picType="png",
                      saveResult=TRUE)

GO Tree at "../results/dge/age-effect/age_type_gender/go/104w_VS_52w_M_GF.png"

104w VS 8w (F, GF)

In [47]:
col = get_GO_network_col(BP_network, "104w VS 8w (F, GF)")
dotRes = getAmigoTree(goIDs=names(col),
                      color=col,
                      filename="../results/dge/age-effect/age_type_gender/go/104w_VS_8w_F_GF",
                      picType="png",
                      saveResult=TRUE)

GO Tree at "../results/dge/age-effect/age_type_gender/go/104w_VS_52w_F_GF.png"

Cellular components

Dot-plot with the most over-represented CC GO (20 most significant p-values for the different comparison)

In [48]:
plot_top_go(age_type_gender_deg,
            "CC",
            40)

Using term, id as id variables
Using term, id as id variables

Molecular functions

Dot-plot with the most over-represented MF GO (20 most significant p-values for the different comparison)

In [49]:
plot_top_go(age_type_gender_deg,
            "MF",
            40)

Using term, id as id variables
Using term, id as id variables

KEGG pathways

In [ ]:
plot_kegg_pathways(age_type_gender_deg$over_represented_KEGG[,"category"],
                   age_type_gender_deg$fc_deg,
                   "../results/dge/age-effect/age_type_gender/kegg/over_repr_kegg/")
In [ ]:
plot_kegg_pathways(age_type_gender_deg$under_represented_KEGG[,"category"],
                   age_type_gender_deg$fc_deg,
                   "../results/dge/age-effect/age_type_gender/kegg/under_repr_kegg/")

What happen to SPF aging genes in GF?

Genes:

  1. Set 52w vs 8w: 52w != 8w (SPF)
  2. Set 104w vs 52w: 104w != 52w (SPF)
  3. Set all: 52w != 8w (SPF) and 104w != 52w (SPF)

Female

Genes in "../results/dge/age-effect/age_type_gender/spf_f_aging_genes"

In [50]:
spf_f_aging = list()
spf_f_aging$deg_52w_vs_8w = rownames(age_type_gender_deg$fc_deg[!is.na(age_type_gender_deg$fc_deg[,"52w VS 8w (F, SPF)"]),])
spf_f_aging$deg_104w_vs_52w = rownames(age_type_gender_deg$fc_deg[!is.na(age_type_gender_deg$fc_deg[,"104w VS 52w (F, SPF)"]),])
spf_f_aging$deg_52w_vs_8w_and_104w_vs_52w = rownames(age_type_gender_deg$fc_deg[!is.na(age_type_gender_deg$fc_deg[,"52w VS 8w (F, SPF)"]) & !is.na(age_type_gender_deg$fc_deg[,"104w VS 52w (F, SPF)"]),])
sapply(spf_f_aging, length)
capture.output(spf_f_aging, file = "../results/dge/age-effect/age_type_gender/spf_f_aging_genes")
deg_52w_vs_8w
1604
deg_104w_vs_52w
280
deg_52w_vs_8w_and_104w_vs_52w
72
In [51]:
col_1 = rep("52w != 8w (SPF)", length(spf_f_aging$deg_52w_vs_8w))
col_1[spf_f_aging$deg_52w_vs_8w %in% spf_f_aging$deg_52w_vs_8w_and_104w_vs_52w] = "52w != 8w (SPF) and 104w != 52w (SPF)"
log2fc_deg_52w_vs_8w = data.frame(SPF_52w_vs_8w = age_type_gender_deg$fc_deg[spf_f_aging$deg_52w_vs_8w, "52w VS 8w (F, SPF)"],
                                  GF_52w_vs_8w = age_type_gender_deg$fc_deg[spf_f_aging$deg_52w_vs_8w, "52w VS 8w (F, GF)"],
                                  GF_104w_vs_52w = age_type_gender_deg$fc_deg[spf_f_aging$deg_52w_vs_8w, "104w VS 52w (F, GF)"],
                                  genes = spf_f_aging$deg_52w_vs_8w,
                                  provenance = col_1)

col_2 = rep("104w != 52w (SPF)", length(spf_f_aging$deg_104w_vs_52w))
col_2[spf_f_aging$deg_104w_vs_52w %in% spf_f_aging$deg_52w_vs_8w_and_104w_vs_52w] = "52w != 8w (SPF) and 104w != 52w (SPF)"
log2fc_deg_104w_vs_52w = data.frame(SPF_104w_vs_52w = age_type_gender_deg$fc_deg[spf_f_aging$deg_104w_vs_52w, "104w VS 52w (F, SPF)"],
                                    GF_52w_vs_8w = age_type_gender_deg$fc_deg[spf_f_aging$deg_104w_vs_52w, "52w VS 8w (F, GF)"],
                                    GF_104w_vs_52w = age_type_gender_deg$fc_deg[spf_f_aging$deg_104w_vs_52w, "104w VS 52w (F, GF)"],
                                    genes = spf_f_aging$deg_104w_vs_52w,
                                    provenance = col_2)
In [52]:
pal = c(rgb(0.5,0.5,1), rgb(0.5,1,0.5,alpha=0.5), rgb(1,0.5,0.5,alpha=0.5))
pal = setNames(pal, c("52w != 8w (SPF)", "104w != 52w (SPF)", "52w != 8w (SPF) and 104w != 52w (SPF)"))

p1 = plot_ly(log2fc_deg_52w_vs_8w,
        x = ~SPF_52w_vs_8w,
        y = ~GF_52w_vs_8w,
        text = paste("Gene: ", log2fc_deg_52w_vs_8w$genes),
        mode = "markers",
        color = ~provenance,
        colors = pal,
        legendgroup = ~provenance) %>%
     layout(xaxis = list(title = "Log2FC 52w vs 8w (SPF)"), yaxis = list(title = "Log2FC 52w vs 8w (GF)"))

p2 = plot_ly(log2fc_deg_104w_vs_52w,
        x = ~SPF_104w_vs_52w,
        y = ~GF_52w_vs_8w,
        text = paste("Gene: ", log2fc_deg_104w_vs_52w$genes),
        mode = "markers",
        color = ~provenance,
        colors = pal,
        legendgroup = ~provenance) %>%
     layout(xaxis = list(title = "104w vs 52w (SPF)"), yaxis = list(title = "52w vs 8w (GF)"))

p3 = plot_ly(log2fc_deg_52w_vs_8w,
        x = ~SPF_52w_vs_8w,
        y = ~GF_104w_vs_52w,
        text = paste("Gene: ", log2fc_deg_52w_vs_8w$genes),
        mode = "markers",
        color = ~provenance,
        colors = pal,
        legendgroup = ~provenance) %>%
     layout(xaxis = list(title = "52w vs 8w (SPF)"), yaxis = list(title = "104w vs 52w (GF)"))

p4 = plot_ly(log2fc_deg_104w_vs_52w,
        x = ~SPF_104w_vs_52w,
        y = ~GF_104w_vs_52w,
        text = paste("Gene: ", log2fc_deg_104w_vs_52w$genes),
        mode = "markers",
        color = ~provenance,
        colors = pal,
        legendgroup = ~provenance) %>%
     layout(xaxis = list(title = "104w vs 52w (SPF)"), yaxis = list(title = "104w vs 52w (GF)"))

subplot(p1, p2, p3, p4, nrows = 2, shareX = TRUE, shareY = TRUE) %>%
  layout(title = "Log2FC") %>% embed_notebook

No trace type specified:
  Based on info supplied, a 'scatter' trace seems appropriate.
  Read more about this trace type -> https://plot.ly/r/reference/#scatter
Warning message:
“Ignoring 1006 observations”No trace type specified:
  Based on info supplied, a 'scatter' trace seems appropriate.
  Read more about this trace type -> https://plot.ly/r/reference/#scatter
Warning message:
“Ignoring 242 observations”No trace type specified:
  Based on info supplied, a 'scatter' trace seems appropriate.
  Read more about this trace type -> https://plot.ly/r/reference/#scatter
Warning message:
“Ignoring 1519 observations”No trace type specified:
  Based on info supplied, a 'scatter' trace seems appropriate.
  Read more about this trace type -> https://plot.ly/r/reference/#scatter
Warning message:
“Ignoring 233 observations”

Male

Genes in "../results/dge/age-effect/age_type_gender/spf_m_aging_genes"

In [53]:
spf_m_aging = list()
spf_m_aging$deg_52w_vs_8w = rownames(age_type_gender_deg$fc_deg[!is.na(age_type_gender_deg$fc_deg[,"52w VS 8w (M, SPF)"]),])
spf_m_aging$deg_104w_vs_52w = rownames(age_type_gender_deg$fc_deg[!is.na(age_type_gender_deg$fc_deg[,"104w VS 52w (M, SPF)"]),])
spf_m_aging$deg_52w_vs_8w_and_104w_vs_52w = rownames(age_type_gender_deg$fc_deg[!is.na(age_type_gender_deg$fc_deg[,"52w VS 8w (M, SPF)"]) & !is.na(age_type_gender_deg$fc_deg[,"104w VS 52w (M, SPF)"]),])
sapply(spf_m_aging, length)
capture.output(spf_m_aging, file = "../results/dge/age-effect/age_type_gender/spf_m_aging_genes")
deg_52w_vs_8w
294
deg_104w_vs_52w
3295
deg_52w_vs_8w_and_104w_vs_52w
119
In [54]:
col_1 = rep("52w != 8w (SPF)", length(spf_m_aging$deg_52w_vs_8w))
col_1[spf_m_aging$deg_52w_vs_8w %in% spf_m_aging$deg_52w_vs_8w_and_104w_vs_52w] = "52w != 8w (SPF) and 104w != 52w (SPF)"
log2fc_deg_52w_vs_8w = data.frame(SPF_52w_vs_8w = age_type_gender_deg$fc_deg[spf_m_aging$deg_52w_vs_8w, "52w VS 8w (M, SPF)"],
                                  GF_52w_vs_8w = age_type_gender_deg$fc_deg[spf_m_aging$deg_52w_vs_8w, "52w VS 8w (M, GF)"],
                                  GF_104w_vs_52w = age_type_gender_deg$fc_deg[spf_m_aging$deg_52w_vs_8w, "104w VS 52w (M, GF)"],
                                  genes = spf_m_aging$deg_52w_vs_8w,
                                  provenance = col_1)

col_2 = rep("104w != 52w (SPF)", length(spf_m_aging$deg_104w_vs_52w))
col_2[spf_m_aging$deg_104w_vs_52w %in% spf_m_aging$deg_52w_vs_8w_and_104w_vs_52w] = "52w != 8w (SPF) and 104w != 52w (SPF)"
log2fc_deg_104w_vs_52w = data.frame(SPF_104w_vs_52w = age_type_gender_deg$fc_deg[spf_m_aging$deg_104w_vs_52w, "104w VS 52w (M, SPF)"],
                                    GF_52w_vs_8w = age_type_gender_deg$fc_deg[spf_m_aging$deg_104w_vs_52w, "52w VS 8w (M, GF)"],
                                    GF_104w_vs_52w = age_type_gender_deg$fc_deg[spf_m_aging$deg_104w_vs_52w, "104w VS 52w (M, GF)"],
                                    genes = spf_m_aging$deg_104w_vs_52w,
                                    provenance = col_2)
#log2fc_deg_52w_vs_8w_and_104w_vs_52w = data.frame(
#    SPF_52w_vs_8w = age_type_gender_deg$fc_deg[spf_m_aging$deg_52w_vs_8w_and_104w_vs_52w, "52w VS 8w (M, SPF)"],
#    SPF_104w_vs_52w = age_type_gender_deg$fc_deg[spf_m_aging$deg_52w_vs_8w_and_104w_vs_52w, "104w VS 52w (M, SPF)"],
#    GF_52w_vs_8w = age_type_gender_deg$fc_deg[spf_m_aging$deg_52w_vs_8w_and_104w_vs_52w, "52w VS 8w (M, GF)"],
#    GF_104w_vs_52w = age_type_gender_deg$fc_deg[spf_m_aging$deg_52w_vs_8w_and_104w_vs_52w, "104w VS 52w (M, GF)"],
#    genes = spf_m_aging$deg_52w_vs_8w_and_104w_vs_52w)
In [55]:
pal = c(rgb(0.5,0.5,1), rgb(0.5,1,0.5,alpha=0.5), rgb(1,0.5,0.5,alpha=0.5))
pal = setNames(pal, c("52w != 8w (SPF)", "104w != 52w (SPF)", "52w != 8w (SPF) and 104w != 52w (SPF)"))

p1 = plot_ly(log2fc_deg_52w_vs_8w,
        x = ~SPF_52w_vs_8w,
        y = ~GF_52w_vs_8w,
        text = paste("Gene: ", log2fc_deg_52w_vs_8w$genes),
        mode = "markers",
        color = ~provenance,
        colors = pal,
        legendgroup = ~provenance) %>%
     layout(xaxis = list(title = "52w vs 8w (SPF)"), yaxis = list(title = "52w vs 8w (GF)"))

p2 = plot_ly(log2fc_deg_104w_vs_52w,
        x = ~SPF_104w_vs_52w,
        y = ~GF_52w_vs_8w,
        text = paste("Gene: ", log2fc_deg_104w_vs_52w$genes),
        mode = "markers",
        color = ~provenance,
        colors = pal,
        legendgroup = ~provenance) %>%
     layout(xaxis = list(title = "104w vs 52w (SPF)"), yaxis = list(title = "52w vs 8w (GF)"))

p3 = plot_ly(log2fc_deg_52w_vs_8w,
        x = ~SPF_52w_vs_8w,
        y = ~GF_104w_vs_52w,
        text = paste("Gene: ", log2fc_deg_52w_vs_8w$genes),
        mode = "markers",
        color = ~provenance,
        colors = pal,
        legendgroup = ~provenance) %>%
     layout(xaxis = list(title = "52w vs 8w (SPF)"), yaxis = list(title = "104w vs 52w (GF)"))

p4 = plot_ly(log2fc_deg_104w_vs_52w,
        x = ~SPF_104w_vs_52w,
        y = ~GF_104w_vs_52w,
        text = paste("Gene: ", log2fc_deg_104w_vs_52w$genes),
        mode = "markers",
        color = ~provenance,
        colors = pal,
        legendgroup = ~provenance) %>%
     layout(xaxis = list(title = "104w vs 52w (SPF)"), yaxis = list(title = "104w vs 52w (GF)"))

subplot(p1, p2, p3, p4, nrows = 2, shareX = TRUE, shareY = TRUE) %>%
  layout(title = "Log2FC") %>% embed_notebook

No trace type specified:
  Based on info supplied, a 'scatter' trace seems appropriate.
  Read more about this trace type -> https://plot.ly/r/reference/#scatter
Warning message:
“Ignoring 191 observations”No trace type specified:
  Based on info supplied, a 'scatter' trace seems appropriate.
  Read more about this trace type -> https://plot.ly/r/reference/#scatter
Warning message:
“Ignoring 3215 observations”No trace type specified:
  Based on info supplied, a 'scatter' trace seems appropriate.
  Read more about this trace type -> https://plot.ly/r/reference/#scatter
Warning message:
“Ignoring 197 observations”No trace type specified:
  Based on info supplied, a 'scatter' trace seems appropriate.
  Read more about this trace type -> https://plot.ly/r/reference/#scatter
Warning message:
“Ignoring 1785 observations”

DEGs 104w vs 8w

SPF, F

In [12]:
comp = "104w VS 8w (F, SPF)"
F_SPF_104w_8w_DEG = rownames(age_type_gender_deg$sign_fc_deg)[!is.na(age_type_gender_deg$sign_fc_deg[,comp])]
F_SPF_104w_8w_DE_FC = age_type_gender_deg$sign_fc_deg[F_SPF_104w_8w_DEG,comp]
names(F_SPF_104w_8w_DE_FC) = F_SPF_104w_8w_DEG
# order by FC
F_SPF_104w_8w_DEG = names(sort(F_SPF_104w_8w_DE_FC))
# selection col to plot
col_to_plot = list("SPF 8w F" = grep("SPF_8w_F_+", colnames(norm_counts), perl=TRUE, value=TRUE),
                   "SPF 104w F" = grep("SPF_104w_F_+", colnames(norm_counts), perl=TRUE, value=TRUE),
                   "GF 104w F" = grep("GF_104w_F_+", colnames(norm_counts), perl=TRUE, value=TRUE),
                   "GF 104w M" = grep("GF_104w_M_+", colnames(norm_counts), perl=TRUE, value=TRUE))
main_title = "DE genes in 104w vs 8w (SPF, F)"

Heatmap of the Z-scores of the DEGs (ordered by log2 FC of the comparison 104w VS 8w (F, SPF))

In [13]:
# get DEG data
deg_data = z_scores[F_SPF_104w_8w_DEG, unlist(col_to_plot)]
# plot heatmap
pheatmap(deg_data,
         cluster_rows=F,
         cluster_cols=F,
         show_rownames=F,
         show_colnames=F,
         annotation_col=gta_annot_col,
         annotation_row=NULL,
         annotation_colors = NULL,
         color=rev(brewer.pal(11, "RdBu")),
         breaks = seq(-3.5, 3.5, length=11),
         main = main_title)

Heatmap of the Z-scores with hierarchical clustering

In [14]:
# get non DEG data and do hierarchical clustering on them
hc = hclust(dist(deg_data), method = "complete")
c_deg_data = deg_data[hc$order,]
# plot heatmap
pheatmap(c_deg_data,
         cluster_rows=F,
         cluster_cols=F,
         show_rownames=F,
         show_colnames=F,
         annotation_col=gta_annot_col,
         annotation_row=NULL,
         annotation_colors = NULL,
         color=rev(brewer.pal(11, "RdBu")),
         breaks = seq(-3.5, 3.5, length=11),
         main = main_title)
In [15]:
mean_z_scores = cbind(apply(z_scores[F_SPF_104w_8w_DEG,col_to_plot[[1]]],1,mean),
                      apply(z_scores[F_SPF_104w_8w_DEG,col_to_plot[[2]]],1,mean),
                      apply(z_scores[F_SPF_104w_8w_DEG,col_to_plot[[3]]],1,mean),
                      apply(z_scores[F_SPF_104w_8w_DEG,col_to_plot[[4]]],1,mean))
colnames(mean_z_scores) = names(col_to_plot)
head(mean_z_scores)
SPF 8w FSPF 104w FGF 104w FGF 104w M
Ptprh-0.3965247-0.4559471-0.4559471-0.4559471
Gprc5c-0.3654629-0.3654629-0.3654629-0.3484024
Rho 2.6036670-0.3926242-0.3926242-0.3926242
Gnat1 2.7597362-0.3875829-0.3875829-0.3851476
Rs1 2.0902734-0.4163774-0.4163774-0.4163774
Rbp3 2.3786827-0.3913755-0.3913755-0.3913755

Heatmap of the mean Z-scores per groups of the DEGs (ordered by log2 FC of the comparison 104w VS 8w (F, SPF))

In [16]:
col_annot = data.frame(age = c("8w",rep("104w",3)),
                       type = c(rep("SPF",2),rep("GF",2)),
                       gender = c(rep("F",3),"M"))
rownames(col_annot) = c("SPF 8w F","SPF 104w F","GF 104w F","GF 104w M")
pheatmap(mean_z_scores,
         cluster_rows=F,
         cluster_cols=F,
         show_rownames=F,
         show_colnames=F,
         annotation_col=col_annot,
         annotation_row=NULL,
         annotation_colors = NULL,
         color=rev(brewer.pal(11, "RdBu")),
         breaks = seq(-3.5, 3.5, length=11),
         main = main_title)

Heatmap of the mean Z-scores per groups of the DEGs (hierarchical clustering)

In [18]:
# get non DEG data and do hierarchical clustering on them
hc = hclust(dist(mean_z_scores), method = "complete")
c_mean_z_scores = mean_z_scores[hc$order,]
#
pheatmap(c_mean_z_scores,
         cluster_rows=F,
         cluster_cols=F,
         show_rownames=F,
         show_colnames=F,
         annotation_col=col_annot,
         annotation_row=NULL,
         annotation_colors = NULL,
         color=rev(brewer.pal(11, "RdBu")),
         breaks = seq(-3.5, 3.5, length=11),
         main = main_title)

Plot of the mean Z-score for GF 104w F and GF 104w M samples in function of the log2 FC of the comparison 104w VS 8w (F, SPF)

In [113]:
data = mean_z_scores
x = F_SPF_104w_8w_DE_FC[F_SPF_104w_8w_DEG]
plot(c(min(x), max(x)), c(min(data), max(data)), type= "n", xlab = "Log2 FC of DEG in 104w vs 8w (SPF, F)", ylab = "Mean Z-score")
points(x, data[, 1], col = 1, pch = 20)
points(x, data[, 2], col = 2, pch = 20)
points(x, data[, 3], col = 3, pch = 20)
points(x, data[, 4], col = 4, pch = 20)
legend("topright", col = 1:4, colnames(data), pch = 19)

SPF, M

In [19]:
comp = "104w VS 8w (M, SPF)"
M_SPF_104w_8w_DEG = rownames(age_type_gender_deg$sign_fc_deg)[!is.na(age_type_gender_deg$sign_fc_deg[,comp])]
M_SPF_104w_8w_DE_FC = age_type_gender_deg$sign_fc_deg[M_SPF_104w_8w_DEG,comp]
names(M_SPF_104w_8w_DE_FC) = M_SPF_104w_8w_DEG
# order by FC
M_SPF_104w_8w_DEG = names(sort(M_SPF_104w_8w_DE_FC))
# selection col to plot
col_to_plot = list("SPF 8w M" = grep("SPF_8w_M_+", colnames(norm_counts), perl=TRUE, value=TRUE),
                   "SPF 104w M" = grep("SPF_104w_M_+", colnames(norm_counts), perl=TRUE, value=TRUE),
                   "GF 104w M" = grep("GF_104w_M_+", colnames(norm_counts), perl=TRUE, value=TRUE),
                   "GF 104w F" = grep("GF_104w_F_+", colnames(norm_counts), perl=TRUE, value=TRUE))
main_title = "DE genes in 104w vs 8w (SPF, M)"

Heatmap of the Z-scores of the DEGs (ordered by log2 FC of the comparison 104w VS 8w (M, SPF))

In [20]:
# get DEG data
deg_data = z_scores[M_SPF_104w_8w_DEG, unlist(col_to_plot)]
# plot heatmap
pheatmap(deg_data,
         cluster_rows=F,
         cluster_cols=F,
         show_rownames=F,
         show_colnames=F,
         annotation_col=gta_annot_col,
         annotation_row=NULL,
         annotation_colors = NULL,
         color=rev(brewer.pal(11, "RdBu")),
         breaks = seq(-3.5, 3.5, length=11),
         main = main_title)

Heatmap of the Z-scores with hierarchical clustering

In [21]:
# get non DEG data and do hierarchical clustering on them
hc = hclust(dist(deg_data), method = "complete")
c_deg_data = deg_data[hc$order,]
# plot heatmap
pheatmap(c_deg_data,
         cluster_rows=F,
         cluster_cols=F,
         show_rownames=F,
         show_colnames=F,
         annotation_col=gta_annot_col,
         annotation_row=NULL,
         annotation_colors = NULL,
         color=rev(brewer.pal(11, "RdBu")),
         breaks = seq(-3.5, 3.5, length=11),
         main = main_title)
In [22]:
mean_z_scores = cbind(apply(z_scores[M_SPF_104w_8w_DEG,col_to_plot[[1]]],1,mean),
                      apply(z_scores[M_SPF_104w_8w_DEG,col_to_plot[[2]]],1,mean),
                      apply(z_scores[M_SPF_104w_8w_DEG,col_to_plot[[3]]],1,mean),
                      apply(z_scores[M_SPF_104w_8w_DEG,col_to_plot[[4]]],1,mean))
colnames(mean_z_scores) = names(col_to_plot)
head(mean_z_scores)
SPF 8w MSPF 104w MGF 104w MGF 104w F
Hist1h2ad1.5357470 -0.3671071-0.3671071-0.3671071
Tfec0.5897554 -0.4613623-0.4635082-0.4705031
Scel1.0443712 -0.4081939-0.4205369 0.3031379
Hist1h3a1.4941967 -0.4527986-0.4876769-0.4876769
Hist1h3e2.5841942 -0.4641478-0.5026988-0.5127714
Hist1h1d1.8757543 -0.4664067-0.4766825-0.4699773

Heatmap of the mean Z-scores per groups of the DEGs (ordered by log2 FC of the comparison 104w VS 8w (F, SPF))

In [23]:
col_annot = data.frame(age = c("8w",rep("104w",3)),
                       type = c(rep("SPF",2),rep("GF",2)),
                       gender = c(rep("M",3),"F"))
rownames(col_annot) = c("SPF 8w M","SPF 104w M","GF 104w M","GF 104w F")
pheatmap(mean_z_scores,
         cluster_rows=F,
         cluster_cols=F,
         show_rownames=F,
         show_colnames=F,
         annotation_col=col_annot,
         annotation_row=NULL,
         annotation_colors = NULL,
         color=rev(brewer.pal(11, "RdBu")),
         breaks = seq(-3.5, 3.5, length=11),
         main = main_title)

Heatmap of the mean Z-scores per groups of the DEGs (hierarchical clustering)

In [24]:
# get non DEG data and do hierarchical clustering on them
hc = hclust(dist(mean_z_scores), method = "complete")
c_mean_z_scores = mean_z_scores[hc$order,]
#
pheatmap(mean_z_scores,
         cluster_rows=F,
         cluster_cols=F,
         show_rownames=F,
         show_colnames=F,
         annotation_col=col_annot,
         annotation_row=NULL,
         annotation_colors = NULL,
         color=rev(brewer.pal(11, "RdBu")),
         breaks = seq(-3.5, 3.5, length=11),
         main = main_title)